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In this paper, we address the problem of voltage regulation in power distribution networks with 
deep-penetration of distributed energy resources (DERs), e.g., renewable-based generation, and storage- 
(-h ■ capable loads such as plug-in hybrid electric vehicles. We cast the problem as an optimization program, 

where the objective is to minimize the losses in the network subject to constraints on bus voltage 
magnitudes, limits on active and reactive power injections, transmission line thermal limits and losses. 
We provide sufficient conditions under which the optimization problem can be solved via its convex 
relaxation. Using data from existing networks, we show that the conditions are expected to be satisfied 
£N) ' by most networks. We also provide an efficient distributed algorithm to solve the problem. The algorithm 

<n : 

ly-j i is asynchronous, with a communication topology that is the same as the electrical network topology. 

We illustrate the algorithm's performance in the IEEE 34-bus and the 123-bus feeder test systems. 

o 

(N , 

I. Introduction 

Electric power distribution systems will undergo radical transformations in structure and 
functionality due to the advent of initiatives like the US DOE Smart Grid, and its European 
counterpart Electricity Networks of the Future. These transformations are enabled by the in- 
tegration of i) advanced communication and control, ii) renewable-based variable generation 
resources, e.g., photovoltaics (PVs), and iii) new storage-capable loads, e.g., plug-in hybrid 
electric vehicles (PHEVs); these distributed generation and storage resources are commonly 
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referred to as distributed energy resources (DERs). It has been acknowledge (see, e.g., [QQ|) 
that massive penetration of DERs in distribution networks is likely to cause voltage regulation 
problems due to the fact that typical ratios of line reactance to resistance are such that bus 
voltage magnitudes are fairly sensitive to variations in active power injections (see, e.g., [|2)). 
Similarly massive penetration of PHEVs can potentially create substantial voltage drops [3]. The 
objective of this paper is to address the problem of voltage regulation in distribution networks 
with deep penetration of DERs. 

As of today, voltage regulation in distribution networks is accomplished through tap-changing 
under- load transformers, set voltage regulators, and fixed /switched capacitors. However, these 
equipments are not designed to deal with rapid voltage fluctuations induced by variations in 
active power generated and/or consumed by DERs. An alternative to the use of these conventional 
voltage regulation devices is to utilize the power electronics interface of the DERs themselves 
to provide reactive power for voltage control, e.g., power electronics-interfaced PV solar rooftop 
installations in residential and commercial buildings. While the primary function of the power 
electronics-based interface is to control active power injections, when properly controlled, they 
can also provide reactive power; thus, they provide a mechanism to control reactive power 
injections, which in turn can be used for voltage control (see, e.g., 01, 0). In this regard, there 
are already existing PV rooftop and pole-mount solutions that provide such functionality (see, 
£-g-> [HI, 0); these solutions are endowed with wireless [6], and power-line communication 
devices Q> which is a key for controlling a large number of devices without having to overlay 
a separate communication network. Additionally, as discussed earlier, the voltage profile of a 
distribution network is sensitivity to changes in active power injections; thus, storage-capable 
DERs, e.g., PHEVs, and demand response resources (DRRs), provide a second voltage control 
mechanism as they can be used, to some extent, to shape active power injections. 

In this paper, we propose a method to solve the voltage regulation problem in distribution 
networks that relies on the utilization of reactive-power capable DERs, and to some extent, on 
the control of active power injections enabled by storage-capable DERs and DRRs. We cast 
the problem as an optimization program where the objective is to minimize the losses in the 
network subject to i) constraints on bus voltage magnitudes, ii) upper and lower limits on active 
and reactive power injections, iii) upper limits on transmission line flows, and iv) upper limits on 
transmission line losses. The decision variables are the bus voltages; however the actual control 
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mechanism to fix these are active and reactive power injections on each bus. This optimization 
problem can be written as a rank-constrained semidefinite program (SDP [0) as observed in 
[|9l, where the decision variable is a positive semidefinite matrix constrained to have rank 1, 
rendering, in general a non-convex problem. The problem can be convexified by dropping the 
rank-1 constraint; the conundrum is to establish when the solution of the convexified problem 
provides also a global solution to the original non-convex problem. 

Recently there has been a sequence of papers on attempting to answer this question spurred 
by the observation in [fTOll that the convex relaxation observed above is tight for many IEEE 
benchmark networks. Several independent works f[TT|- [[T3l provided a partial answer: the convex 
relaxation is tight if the network has a tree topology and certain constraints on the bus power 
injections are satisfied. These results were unified and strengthened in lfl4ll through a investigation 
of the underlying geometry of the optimization problem. 

The first contribution of this paper is in extending the results in 031 to include reactive powers, 
providing i) sufficient conditions under which the optimal solution of the convexified problem 
is also the solution of the original problem, and ii) a simple method to establish whether or not 
the original problem has a feasible solution based on the solution of the convexified problem. 
However, it is important to note that even if the convexified SDP provides the solution to the 
original optimization problem, existing algorithms for solving SDPs are not computationally 
efficient for solving large problems [0; therefore these algorithms are not practical for realizing 
our ideas in a realistic power distribution network with thousands of buses. Furthermore even if 
there is a centralized solver with sufficient computing power, the communication infrastructure 
may not be able to reliably transmit all problem data to one particular location with small 
enough delays; this is where the second contribution of this paper lies — the development of an 
asynchronous distributed algorithm for efficiently solving convexified SDPs on tree networks. 

There is previous work that has addressed the voltage control problem by casting it as an 
optimization problem. For example [15] proposes the use of reactive-power capable DERs for 
voltage control; the objective is to minimize the DER reactive power contributions subject to 
the power flow equations and other constraints; and the solution approach relies on linearizing 
the power flow equations around some operating point, rendering a linear program. In lfl6l . the 
authors also use reactive-power capable DERs for voltage control. They partition the system 
buses into groups and provide an algorithm to solve the optimization problem by breaking it 
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into subproblems, each of which is associated with a particular group of buses — the solution 
that this algorithm provides is not necessarily the globally optimal one. In ifTTl . a centralized 
decision maker maintains a database of historical observations and the corresponding solutions 
to the optimization problem obtained offline. 

The remainder of this paper is organized as follows. In Section [III we formulate the voltage 
regulation problem as an optimization problem. Section |nl] states the main theoretical result of 
the paper, while Section [IV] provides a sketch of the proof (the complete proof is provided in 
the Appendix). Section [IV] may be skipped without lose of continuity for the reader mainly 
interested in the algorithm. Section [V] proposes a distributed algorithm to solve the voltage 
control problem, the performance of which is illustrated in Section [VO Concluding remarks are 
presented in Section IVIII 

II. Preliminaries and Problem Formulation 

In this section, we first introduce the power network model, including relevant notation used 
throughout the paper, and then formulate the voltage regulation problem. 

A. Power System Distribution Model 

Consider a power distribution network with 1,2, ... ,n buses. As of today, such networks are 
mostly radial with a single source of power injection referred to as the feeder (see, e.g., [18]). 
Thus, the network topology can be described by a connected tree, the edge set of which is 
denoted by 8, where (i, k) E 8 if i is connected to k by a transmission line; we write i ~ k if 
bus i is connected to bus k and i ^ k otherwise, and (i, k) to denote a line connected between 
i and k. 

Let Vi = \Vi\Z9i describe the voltage in bus i and define the corresponding bus voltage vector 
as v = [Vi V 2 ■ ■ ■ V n ] T G C n . Similarly, let Pj and Q { denote, respectively the active and reactive 
power injections in bus i, and define the corresponding active and reactive power injection vectors 
as p = [Pi P 2 _ ■ ■ P n } T , q=[QiQ 2 ■■■ Qn] T - Let y ik = g ik - jb ik , with b ik > g ik > 0, denote 
the admittance^ of line (i, k), and let y i{ = jbu denote bus i shunt admittance. Then, the power 

'We adopt the standard assumption that, in normal operating conditions, lines are inductive and the inductive effects dominate 
resistive effects (see, e.g., 1191 . Additionally, although the convention is to write — gik + jbik, we flip the sign of the 
imaginary part as it simplifies the developments. 
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flow equations can be compactly written as 

p + jq = Re{diag(vv^Y H )} + j'Im{diag(vv^Y H )}, (1) 

where Y = [Y ik ] is the bus admittance matrix (sometimes we use A[i, k] to denote the (i, k) th 
entry of a matrix A), (Y^) denotes the Hermitian transpose of v (Y), and diag(-) returns 
the diagonal of a square matrix as a column vector. The active/reactive power flows through 
each line (i,k) are 

Pit = \Vi\ 2 g ik + \Vi\\V k \[b ik sm(9 ik ) - g ik cos(9 ik )], (2) 
Qik = \Vi\ 2 b ik - \Vi\\V k \[g ik sm(9i k ) + b ik cos(6 ik )], (3) 

where 9 ik ■=e i -Q k . 

B. Voltage Regulation Problem Formulation 

We are interested in keeping each individual bus voltage magnitude \Vi\ at some prescribed 
value V- e ^ imposed by performance requirements. In order to achieve so, we rely on a limited 
ability to control active/reactive power injections — this is a realistic scenario in a distribution 
system with deep penetration of renewable-based generation, PHEVs and DRRs — and cast the 
problem as an optimization program where the objective is to minimize the losses in the network. 

Let L ik (V~i,V k ) = P ik + P ki ; then, the total losses in the network are given by L(v) := 
Si k-u k)es Lik^Yi-, Vk) and the voltage regulation problem can be formulated as 

minimize L(v) (4a) 

V 

subject to \Vi\ = V[ ef ,\/i (4b) 

E4<Pi<Pi,Vi (4c) 

Q, < Qi <Qi,Vi (4d) 

\Pik\ <Pik, Vz-A; (4e) 

Lik(Vi, Vk) < Lik, ~ k, (4f) 

where Pj and Qi are functions of v as defined in (OQ), and P ik is a function of Vi and V k as 
defined in ©. Without loss of generality and to ease the notations in subsequent developments, 
hereafter we assume V^ = 1 p.u. for all i. The constraints in (l4bl capture the voltage regulation 
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goal. The constraints in (|4cb and (|4cfl ) describe the limited ability to control active/reactive power 
injections on each bus i; Pi (Pj), and Q { (Q ), respectively denote the upper (lower) limits on 
the amount of active and reactive power that each bus i can provide. The constraints in (|4eb 
capture line thermal limits, while (l4fl ) imposes loss limits on individual lines, this ensures an 
homogeneous distribution of losses across the system. 

The optimization problem in © is difficult for two reasons: i) it is not convex due to 
the quadratic relationship between bus voltages and powers; ii) depending on the size of the 
network, there could potentially be a large number of variables and constraints. Section [III] and 
Section [IV] address i) by convexifying the problem, while Section |V] addresses ii) by proposing 
a computationally efficient algorithm to solve it. 



In this section, we state the main theoretical result, which is that under certain conditions on 
the angle differences between adjacent buses and the lower bounds on reactive power injections, 
the nonconvex problem in © can be solved exactly by solving its convex SDP relaxation. In 
order to state the exact result, it is convenient to write (HJ in matrix form. 

A. Voltage Regulation Problem Formulation in Matrix From 

Let Ej G R nxn , with En = 1 and all other entries equal to zero, and define Aj = |(Y H E; + 
EjY), and Bj = ^-(Y^E; — E,Y). Then, the active and reactive power injections in bus i are, 
respectively, given by Pi = Tr(AjVv^), and Qi = Tr(BiVv^), where Tr(-) is the trace operator. 
For each (i, k) E 8, define a matrix A^, with its (l,m) th entry given by 



III. Convex Relaxation 



A ik [l,m] = < 



-Re{r 4fc } 

Y lk /2 

Yi/2 



if I — m — i 



if / = i and m = k 



if / = k and m = i 



(5) 







otherwise. 



V 
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the active power flow through the (z, k) line is given by P ik = Ti(A ik w\ H ). Let G ik = A ik +A ki . 
Then, we can rewrite © as 

n 

minimize ^^Tr(AjW H ) (6a) 



i=i 



subject to \Vi\ = l,Vi (6b) 

P 4 <TT(A i ynr H ) <P*,Vz (6c) 

Q 4 <Tr(B,vv H ) < Q 4 ,Vz (6d) 

Tr(G ife vv") <I ifc ,Vz~A; (6e) 

|Tr(AjfcW H )| < ~P ik ,\/i ~ fc. (6f) 

Note the outer product vv H is a positive semidefinite rank-1 n x n matrix. Conversely, given a 
positive semidefinite rank-1 n x n matrix, it is always possible to write it as an outer product 
of a vector and itself. Thus, we can rewrite © as 

n 

minimize ^Tr(AjW) (7a) 

i=i 

subject to |W[i,i]| = l,Vi (7b) 

P t <Tr(AiW) <P i? Vz (7c) 

Q. < Tr(BjW) < Qj, Vz (7d) 

Tr(G tfc W H )<I jfc ,Vz~A; (7e) 

|Tr(A ifc W)| <P ik ,Vi~h (7f) 

rank(W) = 1; (7g) 



the rank-1 constraint in (7gi renders the problem non-convex. 
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B. Convexification 

The problem in © is not convex due to the in constraint ( |7g) ; by removing it, we obtain a 
relaxation that is convex: 

n 

minimize ^^Tr(AjW) (8a) 

i=l 

subject to \W[i,i}\ = l,Vi (8b) 

P^Tr^W) <Pi,Vi (8c) 

Q.<Tr(B i W)<Q iJ Vi (8d) 

Ti(G ik W H ) <L ik ,\/i~k (8e) 

|Tr(A ifc W)| <P ik ,Vi~k (8f) 

The convex relaxation is not always tight since the rank of the solution to © could be greater 
than 1. Thus the solution to © does not always coincide with the solution to ©. However, by 
imposing two conditions that are widely held in practice, the non-convex problem in © can be 
solved exactly by solving © as stated in the following theorem. 

Theorem 1. Define 9f k to be the smallest positive solution to the equation P ik = Pi k , with Pi k 
given in ©/or \Vi\ = \V k \ = 1; and 9f k = cos _1 (l - if |jt < 2 and 8g = oo if > 2. 
Let 6i k = mm(9f k , 9f k ). Suppose 6i k satisfies 

- tan" 1 (b ik /g ik ) < 9 ik < tan" 1 (b ik /g ik ) ; (9) 

and reactive power injection lower bounds satisfy 

Qi<Pi, z = 2,...,n, (10) 

with f3i = Y^k-keeu) — 9ik sin(^jfc) — bi k cos(9 ik ), where C(i) is the set of all neighbors of i and 
9{ k = min(tan~ 1 (f iiL ), 9i k ). Let W* be an optimal solution to the relaxed problem ([8]). Then 

1) IfW* is rank 1, then W* = v*(y*) H for some vector v*. Furthermore, v* is the optimal 
solution to the voltage regulation problem stated in ©. 

2) If rank(W*) > 1, then there is no feasible solution to the voltage regulation problem 
stated in ([4]). 

3) If ® is infeasible, then the the voltage regulation problem stated in © is infeasible. 
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Fig. 1. The active line flow region Tih, the reactive flow region Qik, and the linear transformation between them. 

The theorem is proved for a two-bus network in Section |IV]by studying the geometry of the 
feasibility set of the original problem in © and that of its convex relaxation in ([8]). The intuition 
and geometric insight developed by studying the two-bus network carries over to a general tree 
network and the full proof is provided in the Appendix. 

IV. Sketch of Theorem [T] Proof 

The insights into Theorem 1 are obtained by studying the geometry of the sets that result 
from the constraints on line power flows and power injections as described in (|4cT>— d4fb . This 
geometric view was explored in previous works [fTT|. [fT4|. Here, we revisit the results of [fl4l 
and generalize them to include reactive power. 



A. Active and reactive line flow regions 



First, recall from dTbl that \Vi\ = \V k \ = 1 p.u. Then, let £ M 2 and Q^. £ M 2 , respectively, 
denote the regions that contain all the [P ik , P k j\ T and [Qik, Qki] T that can be achieved from © 
and © by varying 9 ik between and 2tt; it is easy to see that for < 9 ik < 2n, © and © are 
linear transformations of a circle. Thus, as depicted in Fig. [Q the active and reactive line flow 
regions 7 ik and Q ik are two ellipses. The center of 7 ik (Q ik ) is [g ik ,g ik } T ([b ik ,b ik ] T ), its major 
axis is parallel to [1, — 1] T ([1, 1] T ) and has length b ik (b ik ), while its minor axis is parallel to 
[1, 1] T ([1, — 1] T ) and has length g ik (gi k ). Both ellipses are related by a linear invertible mapping: 
Qik = HikJ^ik, with 



H 



ik 



2b ik g ik 



bfk-ik b 2 ik + g 



ik 



b 2 ik + 9fk b 2 ik -g lk 



(11) 
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Fig. 2. The flow region under thermal loss constraints (left) and line flow constraints (right). The bold curves indicate the 
feasible part. 



The line flow constraints in (|4eT ) and the thermal loss constraints in (l4fT) appear as linear 
constraints on the line flow regions as shown in Fig. |2] Thus, for each line (i, k), we can replace 
both constraints by a single one, which has the form of the line flow constraint for properly 
defined upper limits. We adopt this convention in subsequent developments. Furthermore, since 
the ellipses have empty interior, this flow constraint can be translated into angle constraints 
on 9 ik of the form \0 ik \ < 9 ik . Conversely, an angle constraint on 9 ik can be converted into a 
flow constraint. Let JFe,%k an d Ge,ik denote, respectively, line (i, k) angle-constrained active and 
reactive line flow regions, then 

Fe,ik = {[Pik, Pki\ T ■ Pik = gik[± - cos(6y + b ik sm(9 ik ), 
Pki = 9ik[i ~ cos(0 ifc ] - b ik sm(9 ik ), \9 ik \ < 9 ik ], 

Qe,ik = {[Pik,Pki} T ■ Pik = b ik [l - cos(9 ik \ - g ik sm(9 ik ), 
Pki = b ik [l - cos(9 ik ] + b ik sm(9 ik ), \9 ik \ < 9 ik }. 

B. Feasible Region of a Two-Bus Network 

Consider a system with only two buses connected by a line (1, 2), with Pi (Qi), and P 2 (Q2) 
denoting the active (reactive) power injections on bus 1 and 2, respectively; for a two-bus system, 
Pi = P12 (Qi = Q12) and P 2 = P21 (Q2 — Q2i)- The relaxed problem in © convexities the 
feasible region of the problem © by filling up the corresponding ellipses as shown in Fig. [3] 
The condition in © is a constraint on the maximum angle difference across the line. Intuitively, 
this angle constraint is such that only the lower left part of the line flow ellipses is feasible. 
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Fig. 3. Angle-constrained active line flow region and its convex hull (filled in blue region) when relaxation is tight (left) and 
when relaxation does not provide solution to the original problem (right). 




Fig. 4. Active power injection region (left) and reactive power injection region (right) under reactive power injection lower 
bound. 



For example, this condition is satisfied by the angle-constrained regions in Fig. |3] Fig. |3] shows 
the intersection of bus power constraints with the angle-constrained active injection regions. In 



Fig. |3(a)[ both bus power constraints are upper bounds. Since the optimal solution of the power 
loss problem occurs in the lower left corner, the convex relaxation is tight. In Fig. |3(b)[ both 
bus power constraints are lower bounds; in this case the optimal solution is inside the ellipses, 
therefore rank W* = 2 and the original problem is infeasible. 

It is important to note that the observations made in Fig. |3]hold as long as the angle-constrained 
injection region only includes the lower left half of the ellipse (as described by ®). From 
thermal data for some common lines in |[T8l . we expect that the angle to be constrained to 
6 ik E [—10°, 10°]. Even for a relatively small b ik /gik ratio of 2, ik = tan -1 -(b ik / ' gik) = 63.4° 
and the condition \9 ik \ < 6 ik is always satisfied. Therefore in most practical networks, it is 
expected that the thermal constraints in the network is small enough that the condition in © 
should be satisfied almost always. 

The second condition in (flOl) is to ensure that the reactive lower bound is large enough such 
that Qi > Q. for all feasible Qi . If the reactive lower bounds are tight at the optimal solution of 
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the relaxed problem, then the rank of the optimal matrix W* is not necessarily 1. Figured shows 
the reason that the condition on the reactive power lower bounds are needed. Figure |4(b)| gives 
the reactive injection region with a tight reactive lower bound on bus 2. Figure |4(a)| shows the 



corresponding active power injection region. Observe that it is possible for the optimal solution 
of the relaxed problem to be rank 2, while the original problem remains feasible. The condition 
Q i < fa rules out this phenomenon by ensuring that the reactive power lower bounds are never 
tight. 

C. General Tree Networks 

The geometrical intuition developed for the two-bus network carries over to a general tree 
network due to two reasons: the flows on each line are independent (no cycles) and the active and 
reactive power injections can be described, respectively, as linear combinations of the active and 
reactive line flows; these are the main ideas used in the proof of Theorem Q] (see the Appendix). 

V. A Distributed Algorithm for Solving the Convexified Problem 

In Section HTl] we showed that the SDP program in ([8]) is a convex relaxation of the voltage 
support problem in ©. Since the objective is to regulate the voltages in the presence of renewable 
sources, the optimization problem needs to be solved no slower than the time scale at which 
the renewable-based power injections change significantly. General-purpose SDP solvers scale 
poorly as the problem size increases [8]. Thus for large distribution networks with hundreds or 
thousands of buses, solving the SDP problem in a minute to sub-minute scale is challenging. 
Furthermore standard solvers for SDP problems are centralized; i.e., all problem data are assumed 
to be presented at a single processor. However the communication infrastructure in a distribution 
network may not be able to transmit the data from all users to a centralized location fast 
enough. By exploiting the tree structure of the distribution network, we propose a distributed and 
asynchronous algorithm to solve (H) that only requires communication between two neighboring 
buses in the network. An independent but related work has appeared in 



A. Algorithm Derivation 

The proposed algorithm consists of two stages: local optimization and consensus. In the local 
optimization stage, each node solves its own local version of the problem. In the consensus 
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stage, neighboring nodes exchange Lagarangian multipliers obtained from the solutions to their 
corresponding local optimums, with the goal of equalizing the phase angle differences across a 
line from both of its ends. 

Let Mi be the set of buses directly connected to bus i by transmission lines, together with bus 
i itself, i.e., Mi = {k : k ~ i, Wk} U {i}. For a n x n matrix M, let M (i) denote the |jVi| x |jV»| 
submatrix of M whose rows and columns are indexed according to Mi. Similarly, for the n X 1 
vector v, yW is the corresponding A/rdimensional vector indexed by Mi. We can rewrite dH) as 



w 



minimize V Tr(A w W (i) ) 

/■W,...,w(«)^o'f-^ 

i=i 

subject to diag(W w ) = v w o v w , Vz 
£i < Tr(A w W w ) < P;,Vi 
Q. < Tr(B (i) W (i) ) < Q^Vi 
|Tr(A®W«)|<P ifc ,V(i,fc)e£ 



ik 



W$ = W£ , ^k)e£, 



r(k) 



(12a) 

(12b) 
(12c) 
(12d) 
(12e) 
(12f) 
(12g) 



where o is the Hadamard product. It is easy to verify that (112al ). (I12bl) . (I12cl) . (I12dl) . and (I12el) 
are equivalent to (|8a1> . (l8bl . (l8cl) . (l8dl) . and ©, respectively, as A* in d8c]), B; in (|8d]> . and A ik 
in (l8fT) have non-zero elements only at (i,k), (k,i), VA; ~ i. Since all Mi's are maximal 
cliques, W ^ is tantamount to W (i) ^ 0, Vz EH, E21. Constraints (fT2fh and (fT2g) 



are added 



to ensure that all W^'s coordinate to form W; in other words, V(z, k) E S, the fl^'s computed 
from and W( fc ) should be the same. 

Let Xik be the Lagrangian multiplier of (11 2fb for (i, k) and similarly for ( 12g ). By relaxing 
(I12fl) and (|12g[), the augmented objective function is 



n 

E 

i=i 



Tr (A«W«) + [^{Wi 

(i,k)e£ 



(0 
fc 



ir 



+A fc4 (wii ) -^ ) )]^x: Tr ( AWw 



i=l 



where A^ is also Hermitian, and its (z, £;) <h entry is i) = if % 



i(0 



(13) 



4> + ^ 



if i < A;, and iii) A$ = - \f h if % > k. With (fl3l) . problem can be divided into n 



iW 
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separable subproblems and the ith subproblem corresponds to bus i, defined as follows: 

minimize Tr(A (i) W w ) (14a) 

ww^o 

subject to diag(W w ) = v w o v (i) (14b) 

Pi < Tr(A (i) W (i) ) < P t (14c) 

Q. < Tr(B w W w ) < Q (14d) 

|Tr(A«W»)| <P ik , Vk~z. (14e) 

We denote the feasible region described by (114bl > — (114-eb together with ^ of Subproblem 
i by d. Define gi (X ik ) = inf wWeC .{Tr(A»W«)}. The gradient of - 9i at X lk is W$\ which 
is the (i, k)th element of the optimal W' 1 '* of g,i determined by solving the ith subproblem 
(fl4"l) . Similarly, that of — g k at X ik is —W^*. Therefore, the gradient of — + g k ) is then 
W$* - W^ ] *. Let W$[t] and W^[t] be W$* and W$* determined at time t, respectively. 
By gradient ascent, at time t + 1, we update X ik by 

X ik [t + 1] = X ik [t] + a[t]{wg[t) - W^it}), (15) 

where a[t] > and X ik [t] are the step size and A^ at time t, respectively. The value of X ki [t + 1] 
can be directly computed from A^ + 1] as X ki = A^. The Lagrangian multiplier X ik is only 
denned for the line (i, k) and the two buses at the ends of the edge, i.e., buses % and k, are 
required to manipulate X ik . The purpose of (TT3T) is to make and as close to each other 
as possible with the help of Xi k . Eq. (fT31) can be computed either by bus i or by bus k and it 
is independent of all other buses and edges. Whenever both the zth and A;th subproblems have 
been computed and so and W^' have been updated, then X ik can then be updated by using 

CCD). 

The optimization problem comprised of (fLTt . together with all the constraints (I14bl) - (ll4el) . 
imposed on the subproblems, is a dual problem of (fT2T ). When all X ik s are optimal, will be 
equal to W^ k for all (z, k)'s and thus the duality gap is zero. Accordingly, we can construct the 
optimal W* of problem (7gi from the values of the W^'s. Furthermore, for a 5 > sufficiently 
small, and due to the convexity of (fl"2l) and the zero duality gap, Algorithm \T\ is guaranteed to 
converge to the optimum. 
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Algorithm 1 Asynchronous Distributed Algorithm 

Given a n-bus network 

1. while \W}1 } - W}^ | > S for any (i, k) e £ do 
2. for each bus i (in parallel) do 

3. Given Aifc,Vfc ~ i, solve J 14b 

4. Return W$\VJfc 

5. end for 

6. Given and Wy* , update Aifc with dT5l > (in parallel 
and asynchronously) 

7. end while 




(a) Network structure (b) Convergence curve 

Fig. 5. A five-bus example 

B. Feasibility 

When the buses determine their own limits on active and reactive powers independently, 
an infeasible problem might result, i.e., an empty feasible region. When there exists a central 
authority having all the bus power information, we can check the feasibility easily. Otherwise, 
it is necessary for the buses to declare infeasibility. 

One sufficient condition for infeasibility of the the problem is that there exists an infeasible 
subproblem (|14l) for any bus. If any bus finds an infeasible subproblem, it is sufficient to say 
that the whole problem is infeasible. To proceed further, the bus with an infeasible subproblem 
should adjust its own active and reactive power limits so as to make the subproblem feasible. A 

(i) (k) 

necessary and sufficient condition for infeasibility is that W> k and W^ k never match for some 
(i,k) G £ when Algorithm Q] evolves. If this happens on edge (i,k), either bus i or bus k or 
both constitute the infeasibility. 
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TABLE I 

BUS INFORMATION OF THE FIVE-BUS EXAMPLE 



Bus 


P 


P 


Q 


Q 


V 


1 


5.2844 


-5.4692 


5.5798 


-5.7604 


1.2247 


2 


-0.0648 


-0.0988 


0.5298 





1.1509 


3 


-0.0423 


-0.5828 


0.6405 





1.1103 


4 


-0.0334 


-0.5155 


0.2091 





0.9762 


5 


-0.0226 


-0.4329 


0.3798 





1.1400 



C. Numerical Performance Enhancements 

Consider the five-bus network given in Fig. |5(a)| Assuming that all A^'s are updated at the 
end of each iteration, the progress of Algorithm [TJ (the curve without power flow constraints) 
and the target optimal objective value are shown in Fig. |5(b)| At iteration 20, when we sum the 
objective function values of all the subproblems, the sum still has around 20% difference to the 
optimal one. Even for a small network, it may take a long time for the algorithm to converge 
to the global optimal solution. Next, we provide some enhancements that improve the algorithm 
convergence speed. 

1) Power Flow Constraints: Constraint (1 1 4-eb means that the active power can flow in any 
direction on the edge (i, k) as long as its magnitude does not exceed the limit Pq~. Assume that 
the global optimal solution W* exists. Our decomposition allows us to compute W* k separately 
by buses i and k, in which each bus determines its local version of W* k , e.g., W$ for bus i. Then 
CE5]) brings both W$ and W$ towards W* k by just equalizing W$ and w£\ If the feasible 
regions Cj and Ck are smaller, it will be easier for ([151 to reduce the discrepancy between 
and 

The additional assumption we make is that all buses are net consumers of active power except 
the feeder; that is, Pi < for i — 2, 3, . . . , n. It is because even for deep penetration levels 
of renewable-based electricity generation in distribution systems, it is unlikely that the power 
generated by these resources in a particular bus will compensate the load at that bus; let alone 
they will compensate for the power coming from the feeder. We now know that the active power 
must flow from buses % to k along the edge (i, k) with i < k, i.e., P ik > 0. With this observation, 



March 3, 2013 



DRAFT 



17 



we can re-write (|14el ) as 

< P lk = Tr(A£?WW) < ?ifc, (16) 
< Pu = Tr(A^W«) < 0, (17) 

from the perspectives of buses i and fc, respectively. We can actually replace (I Hel l for (i, k) 
of Subproblem i by (fT6l) and similarly (I14el > for (z, fc) of Subproblem fc by (fT7l) . If we apply 
the same logic to all edges connecting to bus i, we can construct a smaller feasible region t{ 
for Subproblem i. For the edge (i,k), the constructions of C, and C fc can help and 
converge to W* k faster. 

With this modification, the progress of the algorithm for the five-bus example is also depicted 
in Fig. |5(b)[ where we can see that the algorithm converges faster. 

2) Feasible Solution Generation: When the algorithm converges, we have that 

Tr(A w W (i) ) = Tr(A w W w ), Vi, (18) 

which holds when all its associated A^'s are optimal; this is equivalent to have both of the 
following held: 

Tr(A w W w ) = Tr(A w W (i) *) P, = P*, Vi, (19) 

Tr(B (i) W (i) ) = Tr(B (i) W w *) <=> Q< = Q*, Vi. (20) 

In other words, Algorithm [T] tries to find the the optimal active and reactive power pair [P*, Q*] T 
for each bus i by manipulating A^'s defined for the corresponding lines. The more lines are 
connected to a bus (i.e., the more Ajfc's it involves), the more difficult is for (fT9l and (l20l to 
hold. The [Pi, Qi] T pair affects the [P k , Q k } T pair through \ ik . Consider the situation where edge 
(i,k) is the only line connected to bus k except for bus i. When [Pk,Qk] T becomes optimal, 
this helps bus i converge in the sense that this reduces the variations of [Pi, Qi] T induced from 
bus k. When Algorithm Q] evolves, the [Pk, Qk] T of leaf buj^l k converges first as a leaf bus has 
only one edge. Then, we have the buses connected to the leaf buses converged. We continue 
this process and finally go up to the feeder. 

For any leaf node k, we have Pk = Pki and Q k = Qu, where bus i is the only bus connected 
to bus k. When the algorithm evolves, we obtain W^* from the solution of the kth subproblem 

2 A leaf bus is a bus with one line attached only. 
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(fl4l) when Tr(A (fc) W (fc) ) and Tr(B (fc) W (fc) ) are equal to P fe * and Q* k , respectively. Once we 
have fixed W( k , we can add the constraint W^ k = W^ k to the zth subproblem for bus i 
by passing a message containing the value of from bus k to bus i. In matrix form, this 

constraint is equivalent to Tr(C^W^) = Re{W^*} and Tr(D«W«) = ha{W^ ] *}, where 
C« = (c£,Z,m G A/i), with c£ = | if i = i and m — k, C® — \ if I — k and m = i, 
and c£! = otherwise; and D« = {Df m ,l,m E M), with L>£ = ±7 if Z = i and m = k, 
Dim = —\j if I — k and m — i, and I?^ = otherwise. In this case, we reduce the n-bus 
network into the (n — l)-bus one by removing bus k. When all other buses with positive active 
power flown from bus i (i.e. {I : I ~ i, I > i}) have been fixed and "removed", bus i becomes a 
leaf bus in the reduced network. This process continues until we find all W* m ,^(l,m) E S. The 
global solution W* can be constructed from those Wj^'s. However, for any bus k, if we fix P k 
and Qk which is not optimal, these errors will make its connecting bus i being fixed afterwards 
result in incorrect P { and Qi, which are not optimal either. To achieve this, we observe Pi[t] 
and Qi[t] for a certain time period and check if their variations are significant. Assume that we 
are at time t, for the active power, we can keep track of the previous T Pj's and the current 
Pi[t], i.e. [Pi[t - T], Pi[t - T + 1], . . . , Pi[t]] T . We can say that P t [t] has been converged if its 
cumulative change is less than a certain threshold 7 (e.g. 10~ 4 ), i.e., 

^ |p.[ t - r + fc ]-P,[t-T + fc + l]| 

h m-T+k}\ <7 ' (21) 

with a similar condition for the reactive power. 

3) Hot Start: The problem needs to be solved repeatedly; when there are changes to the 
active/reactive limits at any bus, we apply Algorithm Q] to the problem again. In each update, we 
usually have small variation between the new P 4 and the previous ones and also for P i . Thus, in 
subsequent instances of the problem, the optimal angle difference across each line usually does 
not vary significantly. Therefore, we can set \ik [0] with the optimal \* k which can be determined 
from the previous optimal W* k . 

VI. Case Studies 

We test the algorithm performance on the 34- and 123-bus test feeder systems [|23l . Assume 
that the system operates under nominal conditions specified by the given dataset, for one day 
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Fig. 6. Irradiance of a particular day 



(1440 minutes) with constant load demands Pj,Vi. We assume that on each bus i, there are an 
energy storage and connected photovoltaics (PV)-based electricity generation resources, which 
can supply active power, Pf v to the bus locally, i.e. their net effect is to reduce the load. If all 
Pf v is consumed locally, then the active power injection at bus % will be Pj = Pi+Pf v < 0. The 
computed optimal P* £ [p, Pi], i — 2, . . . , n, will then be adjusted by controlling the amount 
of power from the PV devices which will be stored at the local storage. Let be the nominal 
reactive power injection at bus i. It turns out the power electronics on modern PV panels can 
supply reactive power in a range that is sufficient to cancel the nominal reactive power [1241 . so 
we assume that the reactive power can be adjusted in the ranges specified by i) Q { £ [0, 1.2Qi], 
if Qi > 0, and ii) Qi = [— 1.2Qx, 0] otherwise. 

We consider irradiance data in Fig. |6] which correspond to a particular day in November 201 1 
collected at the University of Nevada [1251 : the data are sampled every minute. Assume that the 
PV module at bus i can provide up to 20% of its nominal load p. Thus Pf y, s at different buses 
are different and are proportional to their respective p. At different time, Pf v varies according 
to the irradiance data. 

Since there is only radiation between the 377th and 991th minutes, we define the time 
horizon as [377,991]. We execute the algorithm every one-minute cycle, corresponding to the 
PV generation changes in this period. Recall that the algorithm requires inputs of Lagrangian 
multipliers as the starting points. In minute t, where t £ [377, 991], the inputs to the algorithm are 
the Lagrangian multipliers computed by the algorithm at time t — 1. Moreover each Lagrangian 
multiplier is only stored and manipulated by the two buses at the two ends of the corresponding 
line. The initial Lagrangian multipliers, i.e. at t = 377, are computed from the nominal system 
settings. 
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(a) 37-bus system (b) 123-bus system 

Fig. 7. Objective function values of the solutions computed by the distributed algorithm 
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In order to check if the distributed algorithm can obtain the global optimums, we compare 
the objective function values computed by the distributed algorithm and those by the centralized 
solver; the results are plotted in Fig. [7J For the 37-bus system, we can see that the distributed 
algorithm converges to the optimums all the time, whereas for the 123-bus system convergence 
occurs most of the time. In the 37-bus and 123-bus systems, the centralized solver failed to 
solve the system due to convergence issues. Fig. [8] shows the corresponding computational 
times elapsed. Here we only consider the CPU time spent on the SDP solver and assume that 
communication overheads can be neglected. In our simulation, we implement the algorithm 
iteratively; in each iteration, we solve the subproblems sequentially. In Fig. [8l each subfigure 
contains two curves. One (distributed) is to sum the CPU times of the subproblems which need 
the longest CPU time in each iteration. In other words, we only consider the most demanding 
subproblem in each iteration and then sum the CPU times spent on these subproblems in all 
iterations. The average CPU times for the three cases are 1.28 s, 3.33 s, and 19.69 s, respectively, 
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and they are short enough when each cycle lasts one minute. 

VII. Concluding Remarks 

We proposed to solve a convex optimization based method to solve the voltage regulation 
problem in distribution networks. We cast the problem as a loss minimization program. We 
showed that under broad conditions that are likely to be satisfied in practice, the optimization 
problem can be solved via its convex relaxation. We then proposed a asynchronous and distributed 
algorithm that can be implemented in a network with a large number of buses; we demonstrated 
the effectiveness of the algorithm with two case studies. 

Appendix 

Proof of Theorem [TJ- The first and third cases are clear. The interesting case to prove is to 
show that if W* has rank higher than 1, there are no rank 1 W that is feasible. 

The requirement that Q . < fa for i = 2, . . . , n is to ensure that the reactive lower bound is 
in fact never tight for all the nodes in the network. Let h be parent of i and A; be a child of 
i. Since we assume that power always flow from parent to child in the network, and from the 
angle constraint in ®, < 9m < tan _1 (— ). Over this range, Qih > 0. This corresponds to the 

9ik 

intuition that reactive power should flow up the tree to support the voltage. Note the inductive 
line is very lossy in terms of reactive powers, therefore i might receive or supply reactive power 
to k. The Qik is monotonic in 9q. starting at 9 ik = until it reaches is minimum at an angle 

of tan -1 (g). Let 9 ik = min(tan _1 (g), then Qi = Q lh + Y.k±&c{i) Qik > Y.k±ac{i) Qik > 
T,k:kec(i) b ik ~ 9ik sin(^ fc ) - b ik cos(9(9) ik ) = fa. Therefore if Q. < fa, the reactive power lower 
bounds are never tight. 

To finish the proof we need to introduce some notations from |[T4ll . For a n bus network, 
let {9 ik } be the set of angle constraints, one for each line. Then, the angle-constrained active 
power injection region is the set of all active power injection vectors that satisfy the line angle 
constraints, i.e., V e = {p : p = Re{diag(vv^Y H )}, \Vi\ = 1, \9 ik \ < 9 ik }. Let JF e C R 2n ~ 2 
be the Cartesian product region of the n — 1 active line flow regions, then JRq = H^JFe ,ik- Let 
M G ]R nx2ri - 2 be a matrix with the rows indexed by the buses and the columns indexed by the 
2{n — 1) ordered pair of edges, i.e., if i is connected to k, both (i,k) and (k,i) are included; 
thus M[z, (k, I)] = 1 if i = k or i = I, and M[i, (k, I)] = otherwise. M is a generalized edge 
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to bus incidence matrix, and Vq = MT$ i.e., the power injection region is obtained by a linear 
transformation of the product of line flow regions. 

Similarly, Q e is the product region of the n — 1 reactive line flow regions. Then, for all (i, k), 
by stacking the H, fc 's as defined in (fTTT) into a 2{n — 1) x 2(n — 1) block diagonal matrix H, i.e., 
H = di&g({H ik }i^ k ), we obtain a the global transform between T and Q. The angle-constrained 
reactive power injection region Qg is given by Qq = NlQo = M.HTg. 

Since the lower reactive power injections bounds are never tight by construction, we ignore 
them from now on. Let V be the feasible region of the original problem ©, that is, V = {p : 
3v G C n ,Pt = Tr(A lV v^), \V-\ = l,P i <P i < P i7 Tr(B. t vv") < Q t , \6 ik \ < 6 ik ,\/i ~ k}. We 
can equivalently write V as V = Wl(Te^TpP\To) , where Tp is the flow region satisfying the real 
power constraints, that is, Tp = {f G IR 2n ~ 2 : p = Mf, P t < P, < Pi}. Tq is the flow region 
satisfying the reactive power constraints, that is, Vq — {f G IR 2n_2 : q = MHf , Qi < Q^. 
Since Tp and J-q are defined by linear inequalities, they are convex. However, Tq is not. 

Let S be the feasible region of the relaxed problem d8}. It turns out that S = M(convhull(J r e)n 
Tp n Tq), is convex, and contains V. 

Now we need to define the Pareto-front of a set. Let X C IR n , we say x G X is Pareto-optimal 
if G <Y such that y < x with strict inequality in at least one coordinate. The set of Pareto- 
optimal points is called the Pareto-front of X, and labeled O(X). When minimizing a strictly 
increasing function, the optimal is always achieved in the Pareto-front. Therefore to show the 
second statement in the theorem, it suffices to show the following lemma, 

Lemma 2. Suppose V is not empty, then V = 0(S). 

Suppose the lemma is true, then if the optimal solution of the relaxed problem dHJ is of rank 
2, then V must be empty. 

The proof of this lemma is similar to Lemma 4 in [Q3). Let p* G S be the optimal solution 
of the relaxed problem, f* G convhull(J-6») H Tp H Tq its corresponding active flow vector and 
r* = Hf * be the corresponding reactive power flow vector. It suffices to show that if P* > P is 
then (f* k , /jjfj G Te,ik for every k ~ i. Once this fact is established, the rest of the proof is the 
same as the proof of Lemma 4 in [fT4l . Suppose that P* > P_ t , but (f* k , f ki ) £ Te,ik for some k. 
Then there exists e > such that (f* k - e, /*•) G convhull(J* ei i fe ). Let (f ik , f ki ) = (f* k - e, /£). 



March 3, 2013 



DRAFT 



23 



Since 



— e 




h 2 - n 2 

u ik iJik 

.2 , „2 



+ 9ik 



2bik9ik 

Therefore (fik, fki) is a better feasible flow on the line (i,k), which contradicts the optimality 
off*. ■ 
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